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SYNOPSIS 

A spatial water balance model for the Tarrawarra catchment in Australia has been set up in 
order to investigate the effect of changes in dynamic land uses on water balance. Detailed 
Meteorological and runoff data, soil profiles, soil moisture patterns were available for this 
catchment. The data were processed to prepare the input files required to run the WBM 
model. The simulations were conducted using a 6 minute time step. The model parameters 
were set by a careful analysis of the available field data, with subsequent hydraulic 
conductivity parameters. The model performed reasonably well, although there was some 
tendency for the model to overestimate the surface runoff. This showed difficult and 
challenging task to simulate saturation excess runoff. Different calibration trials were 
attempted to achieve a best combination of model parameters for which there would be 
exact match in time series graphs. Finally a set of parameters were concluded based on 
time series graph and eiTor statistics results. 

Effect of introduction of different types of forests in the catchment showed logical results; 
which implied that the WBM model could be used for prediction of differences in water 
balance and suitably developing an urban area in the catchment. This will provide a 
direction to the city planners to conserve water resources by planting a definite amount of 
land for buildings and selection of appropriate types of trees to be planted. 



CONTENTS 


CERTIFICATE i 

SYNOPSIS ii 

ACKNOWLEDGEMENTS iii 

CONTENTS - ’ iv ' 

LIST OF FIGURES viii 

LIST OF TABLES xi 

Chapter 1 - INTRODUCTION 

1.1 General 1 

1 .2 Role of spatial data tools in hydrologic modeling 3 

1.3 Objectives 5 

1 .4 Layout of the report 5 

Chapter 2 - LITERATURE REVIEW 

2.1 General y 

2.2 The major types of Hydrological models 9 

Chapter 3 ~ THE WBM 

3.1 General , 

3.2 Coupling with GIS 



3.3 Input and output of the WBM 

3.4 Model component 


18 


V 


20 


ChaBt er 4 - STUDY AREA AND DATA PROCESSING 


4. 1 General 29 

ij 

4.2 Desci'iption of the catchment 29 

4.3 The Tarrawarra data set 31 

4.4 Input files used in the model 35 

4.5 Generation of input files from available data 39 

4.5.1 Soil classification 39 

4.5.2 Rainfall 41 

4.5.3 Temperature 42 . 

4.5.4 Evapotranspiration 42 

4.5.5 Runoff 42 

4.6 Otput files obtained from the model 44 

4.7 Data processing 45 

4.8 Calibration and Validation data set 43 


Chapters - PARAMETER ESTIMATION 

5.1 General 

5.2 Calibration using 1995 data set 
5.2.1 Trial 1 


5.2.2 Trial 2 


60 



V A 


5.3 Calibration using 1996 data set 63 

5.3.1 Data of 1996 setl 64 

5.3.2 Data of 96 set2 66 

5.4 Final 1996 trial 

5.4.1 Physics involved in soil moisture modelling 72 

5.4.2 Discussion on Surface runoff and soil moisture profile 77 

// 

5.5 Calibration and validation using newly generated data set 81 

5.6 Results and discussion 85 


Chanter 6 - SENSITIVITY ANALYSIS OF DYNAMIC LAND USE 
PARAMETERS 


General 

89, 

Dynamic land use 

89 

Sensitivity of parameters for grass 

92 

6.3.1 Changes in root depth 

92 

6.3.2 Changes in vegetation cover 

93 

6.3.3 Changes in leaf area index 

94 

6.3.4 Changes in Unproductive land 

95 

Different Land uses 

97 

6.4. 1 Results 

98 

6.4. 1.1 Coniferous forest 

98 

6.4. 1.2 Deciduous forest 

100 


6.4. 1.3 Mixed forest 


101 



Vll 


Chapter 7 - SUMMARY AND CONCLUSIONS 


7.1 Summary 

104 

7.2 Conclusions 

105 

7.3 Scope for future work 

106 

REFERENCES 

108 

Appendix A 

115 

Appendix B 

117 



Vlll 


LIST OF FIGURES 


FIGURE CAPTION 

NO. 

PAGE 

NO. 

3.1 

Vertical column considered for a raster element 

16 

3.2 

Flow chart of “IDRISI - WBM - FLIESS” 

18 

3.3 

Cascade model of raster elements 

21 

3.4 

Schematic representation of a raster element 

22 

3.5 

Transpoitation of rain water through Soil columns 

24 

4.1 

Aerial photo of Tarrawarra catchment looking south east 

30 

4.2 

Tarrawarra catchment with contours lines and 
meteorological stations 

33 

4.3 

Digital elevation model of Tarrawarra catchment 

35 

4.4 

Regions with different soil types 

40 

4.5 

Time series graph of available rainfall 

41 

4.6a 

Time series from the available observed data 

43 

4.6b 

Time series from extracted data 

44 

4.7a 

Time series plot for calibration period 

50 

4.7b 

Time series plot for validation period 

50 

5.1 

Time series graph for trial 1 

56 

5.2a 

Time series graph for Kst of 1 .0 

57 

5.2b 

Time series graph for Kst of 3.0 

58 

5.2c 

Time series graph for Kst 5.0 

58 

5. 2d 

Time series graph for Kst 9.0 

59 

5.3a 

Time series plot for Kst 0.5 and trial 2. 

61 

5.3b 

Time series plot for Kst 0.8 and trial 2 

61 

5.3 c 

Time series plot for Kst TO and trial 2 

62 

5.3d 

Time series plot for Kst 3.0 and trial 2 

62 

5.3e 

Time series plot for Kst 9.0 and trial 2 

63 



IX 


5.4 Time series plot for 96setl 64 

5.5 Time series plot corresponding to minimum error statistics 68 

5.6 Best observed match in time series plot for 96 set 2 69 

5.7 Variation of RMSE with increase in Kf values 70 

5.8 Variation in RMSE and NSC with increase in maximum 71 

infiltration values 

5.9 Variation in MBE with increase in Maximum infiltration 71 

values 

5.10 Processes considered in WBM model 74 

5.11 Infiltration and percolation curve 75 

5.12 Time series plot obtained from final trial 7 8 

5.13 Time series plot obtained from literature 78 

5.14 Soil moisture profile obtained from literature for 80 

Tarrawarra catchment 

5.15 Soil moisture profile obtained from model 80 

5.16a Time series graph for calibration period of Low peak 8 1 

event i.e. period 4L 

5. 16b Time series graph for validation period of Low peak 82 

event i.e. period 2L 

5. 16c Time series graph for calibration period of m edium peak 83 

event i.e. period 4M 

5. 16d Time series graph for validation period of Medium peak 84 

event i.e. period 2M 

5 . 1 6e Time series graph for calibration period of High peak 84 

event i.e. period 4H 

5. 16f Time series graph for validation period of High peak ' 85 

event i.e. period 2H 

6 . 1 Effect of changes in root depth 93 

6.2 Effect of changes in BG 94 

6.3 Effect of changes in LAI 95 



X 


6.4 

Partition of a pixel according to land use 

96 

6.5 

Effect of introduction of unproductive land 

97 

6.6 

Cumulative values of different components of rain for 
coniferous forest 

99 

6.7 

Cumulative values of different components of rain for 
deciduous forest 

101 

6.8 

Cumulative values of different components of rain for 
mixed forest 

102 

A-1 

Typical layers found in a soil profile 

117 



XI 


LIST OF TABLES 


TABLE CAPTION PAGE 

No. NO. 

5.1 ' Best combinations achieved from arbitrary selection 55 

5.2a Error statistics of different calibration trials for 1995 87 

5.2b Error statistics of different calibration trials for 1996 88 

5.2c Error statistics of different calibration and validation for 1996 88 

6.1 Haude monthly coefficient 98 

6.2 Seasonal variations for coniferous forest 99 

6.3 Seasonal variation for deciduous forest 100 

6.4 Seasonal variations for mixed forest 102 



Chapter 1 


Introduction 

1 

1.1 General 

Mathematical models of watershed hydrology are fundamental to water resources 
assessment, development and mangement. They are used to analyse the quantity and 
quality of streamflow, reservoir system operations, groundwater development and 
protection, surface water and groundwater conjuctive use management, water distribution 
systems,_ water use and a range of water resources management activities (Wurbs 1998)1 
Watersheds may be modelled by a lumped model using basin average input data and 
producing total basin streamflow. Such a model may produce reasonable results but 
because of the distributed nature of hydrological properties like soil type, slope and land- 
use, the model cannot be expected to accurately represent the watershed conditions. 
Distributed hydrological models like TOPOG, Water balance model (WBM) provide an 
effective simulation tool for exploring basin hydrological processes and predicting the 
effects of change on catchment response. Distributed models preserve the spatial 
variability by dividing the model domain into small, interconnected cells, which makes 
them suitable for distributed hydrologic forecasting to better represent watershed features. 
The conjunctive use of distributed hydrologic models and field or remote observations 
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permit the comparison of spatial patterns in hydrology which ultimately increase our 
insight on the hydrological processes operating in a basin. 

Modeling of the hydrologic system has been the subject of research among hydrologists 
and engineers for a very long time. The transformation of rainfall into runoff is an 
extremely complex, dynamic, and non-linear process which is affected by many factors 
which are often inter-related. The factors affecting the runoff response of a watershed 
subjected to rainfall input include (a) storm characteristics i.e. intensity and duration of 
rainfall event, (b) watershed characteristics i.e. size, shape, slope, and storage 
characteristics of the watershed, percent of the watershed contributing runoff at the outlet 
at various time steps during a rainfall event, etc., (c) geo - morphological characteristics of 
a watershed i.e. topography, land use patterns, vegetation, soil types, etc., that affect the 
infiltration, and (d) climatic characteristics such as temperature, humidity, wind 
characteristics, etc. The influence of these factors and many of their combinations in 
generating runoff is an extremely complex physical process and is not understood clearly 
(Zhang and Govindaraju, 2000). Many of the rainfall-runoff models rely on the fact that 
the dynamic effects of various factors mentioned above are embedded in the rainfall and 
runoff data. 

Historically, Hydrologists and. researchers have employed deterministic model or 
Conceptual model that considers the physics of the underlying process. These models use 
the equations of mass, energy, and momentum to describe the movement of water over 
land surface and through the unsaturated and saturated zones of the earth (Birikundavyi et 



al., 2002). There are several hydrologic models developed by different organisations across 
the globe. There is no common agreement to accept any of the developed models as the 
standard one due to heterogeneous land use and unpredictable climatic conditions for the 
different catchment systems. The use of emerging technologies of Remote sensing and 
Geographical information system (GIS) has proven their usefulness in development df 
efficient hydrologic models. 

1.2 Role of spatial data tools in hydrologic modeling 

Many types of hydrological analysis are limited by a lack of spatial data, since traditionally 
hydrological data are point measurements. Remote sensing data are fundamentally 
different, they incorporate spatial information. The rapidly developing geographical 
information system (GIS) technology is a powerful tool to organize, process and visualize 
spatial data. The GIS technology can be usefully employed for data management, data 
retrieval, manipulations, organisations and the like - that may be entirely suitable for 
development of distributed models. Thus, remote sensing (RS) and GIS can complement 
each other and enable hydrological models to be more physically based and more efficient. 
It is only during the last five years that an integration of RS and GIS with hydrological 
models became technically feasible (Baumgartner et al 1997). GIS has provided a new 
environment to develop distributed hydrological models. These models take into account 
and predict the values of the studied phenomena at any point within the watershed (Mitas 
and Mitasova 1998). Zhang (1997) developed a framework integrating a stochastic space- 
time rainfall model and a distributed hydrologic simulation with GIS. 
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The remote sensing technology has many advantages (Schultz, 1988). It produces areal 
measurements in place of point measurements. All informations are collected and stored at 
one place. It offers high resolution in space and time. Data are available in digital form. 
Data acquisition does not interface with data observation. Data can be gathered for remote 
areas that are otherwise inaccessible. Once the remote sensing networks are installed, data 
measurements are relatively inexpensive. 

Remote sensing data from Satellite are useful in watershed modelling in two ways. First, 
Satellite data can be used to better define soils and land covers over a watershed, which are 
needed to determine infiltration, evapotranspiration and runoff Second, remote sensing 
measures data over space rather than at a point, and can, therefore be used to correct errors 
in input data based on point measurements as is with rainfall, evaporation. Using these 
techniques, it is possible to generate three dimensional watershed model, including 
precipitation, surface runoff, porous media flow, open channel flow, interaction of ground 
water and channel flow, and transport of water to the atmosphere by evaporation and 
transpiration. 

The GIS can efficiently combination the different data such as topography, soil maps, 
rainfall distribution etc. The integrated use of satellite imagery and GIS technologies 
provide spatial and temporal input parameters to the models. Remote sensing permits the 
observation of data over large and inaccessible areas and monitoring of land use changes in 
a river basin. The effect of land use changes on hydrologic cycle can be understood by 
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estimating the sensitive parameters influencing the hydrological processes. Dynamic land 
use parameters can be thought of in terms of root depth, vegetation density, leaf area index, 
percentage of unproductive land i.e. inrban areas, surface roughness, albedo and 
interception storage. 

1.3 Objectives 

The objective of this work can be summarised in four steps: 

1) To understand the importance of conceptual hydrologic models; 

2) Explore how to use the WBM model and incorporate RS data in simulating 
hydrological processes of a catchment; 

3) Estimate the values of parameters of the WBM model; 

4) Prediction of change in water balance due to dynamic land use changes in the 
catchment using calibrated WBM model. 

1.4 Layout of the report 

This report is specifically intended to give in detail description of using hydrologic 
modelling to observe the change in water balance due to dynamic land use. The presents 
chapter gives an overview of hydrologic modelling and role of RS and GIS in hydrologic 
modelling and presents the objective of the thesis. At the beginning, in second chapter, 
literature review has been presented. Third chapter discusses about the WBM model used 
in this study. Fourth chapter describes the catchment under study along with summary of 
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available data; preparation of input files required to run the model and selection of periods 
for calibration and validation of the model. Next chapter deals with estimation of model 
paiameters using trial and error procedure. Chapter 6 discusses the sensitivity of dynamic 
land use parameters and its affect on water balance. The last chapter presents the 
conclusion of the study and scope for future study. ,, 



Chapter 2 


Literature Review 


2.1 General 

Watershed hydrology is the branch of hydrology which deals with the infiltration of 
hydrologic processes at the watershed to determine the watershed response (Singh, 2002). 
The watershed may be as small as a parking lot or as a large as river basin. The 
quantitative description of the land phase of the hydrologic cycle is very complicated due 
to spatial non-uniformity defined by climate, topography, geology, soils, vegetation and 
land use. Watershed modelling is required for varied purposes such as water resources 
assessment like water quality evaluation, development like planning and designing soil 
conservation practices and management like irrigation water management. These models 
are utilized to quantify the impacts of watershed management strategies, linking human 
activities within the watershed to water quantity and quality of the receiving stream or lake 
for environmental and water resource protection. There is abundance of literature available 
in the area of hydrologic modelling using conceptual methods. This chapter reviews some 
of the important hydrologic models relevant to the present study. 

In the early 1960's the Stanford Watershed Model (SWM) was instrumental in introducing 
the civil engineering profession to the concept of continuous hydrologic modeling. 



8 


Crawford and Linsley (1966) developed SWM to describe entire hydrologic cycle. SWM 
has undergone many revisions since its development initially. Some other popular models 
from the computer age are tlie watershed models of Dawdy and O' Donnell (1965); HEC-1 
(Hydrologic Engineering Center 1968); SWMM (Metcalf and Eddy et al. (1971); PRMS 
(Leavesley et al. 1983); NWS river forecast system (Burnash et al. 1973a,b); SSARR 
(Rockwood 1982); Systeme Hydrologique European (SHE) (Abbott et al. 1986a,b); 
TOPOMODEL (Beven and Kirkby 1979); and IHDM (Morris 1980). Most of these models 
use the laws of physics to simulate various components of the hydrologic cycle. 

During the decades of 1970s and 1980s, with the digital revolution, a number of 
mathematical models were developed not only for simulation of watershed hydrology; but 
also for their applications in other areas, such as environmental and ecosystem 
management. 

Singh (2002) has presented an overview of the hydrologic models used in different 
countries. HEC-HMS is considered as the standard model in the private sector in the 
United states for the design of drainage systems, quantifying the effect of land use change 
on flooding etc. The UBC and WATFLOOD models are popular in Canada for hydrologic 
simulation. The RORB and WBN models are commonly employed for flood forecasting, 
drainage design and evaluating the effect of land-use changes in Australia. TOPMODEL 
and SHE are the standard models for hydrologic analysis in many European countries. The 
HBV model is the standard model for flow forecasting in Scandinavian countries. The 
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ARNO, LCS and TOPIKAPI models are popular in Italy. The Tank models are well 
accepted in Japan. The Xinanjiang model is a commonly used model in China. 

Most of the detenninistic hydrologic mpdels need a large amount data for calibration and 
validation purposes and are computationally extensive. The parameters of the model needs 
to be first calibrated before it can used for predicting runoff. Calibration of the model has 
been the subject of research among hydrologists. Normally, rainfall and runoff data for 
representative storms taken from a watershed under consideration are employed for this 
purpose. New technique of Genetic algorithms is becoming popular in finding the optimal 
parameters of a model. It is interesting to understand modelling concepts and 
interdependency between model parameters and error objective function surfaces, 
particularly when the model is used to determine the likely effects of land use changes. In 
this case, it requires to deal with spatial variability and scaling. 

2.2 The major types of Hydrological models 

Durand [2002] reported that most of the models used in hydrology fit into one of the 
following categories. 

i. Black box empirical models (BBMs): The relationship between the input signal and 
the output signal is derived from long series of observations, using statistical time 
series analysis, multivariate regression or computing tools. 
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ii. Lumped conceptual models (LCMs): The catchment is considered as a whole, and 
the system is modelled using a conceptual, “metaphysical” vision, e.g. a cascade of 
reservoirs or a simplified physical system. 

iii. Semi - distributed, physically - based models ( SDMs); The catchment is divided 
into “fimctional” subparts, which can be areas that are homogeneous with respect to 
one or more characteristics (land use, soil type), or structural elements (slope, 
drainage network), or classes of points having the same hydrological behaviour 
(hydrological similarity concept). The fluxes are computed using physical 
equations with same assumptions and solved analytically. 

iv. Distributed, process based numerical models (DPMs): The catchment is described 
using a grid (regular or not, square or not). The fluxes are computed using the 
partial differential equations of the physics of porous media and hydraulics, solved 
by a numerical scheme (finite differences or finite elements). The WBM model 
used in the present study falls under this category. 

Four important hydrologic models are reviewed to understand the structure of the 
hydrologic models in terms of input and output data and the practical applicability of these 
models. 

TOPOG 


TOPOG is a research tool developed by a large group effort of CSIRO (Commonwealth 
scientific and industrial research organisation) Land and water and the Cooperative 
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Research Center for Catchment Hydrology, Australia, in 1986 for wide range of water and 
land management. It is a Hydrologic modelling package, which simulates the movement of 
water through the complex terrain. It explains the behaviour of hydrologic system based 
on physics behind it and hence, deterministic in nature. The model accounts for distributed 
parameters such as soil types, vegetation and climate by considering their spatial 
variability. It describes the movement of water over the surface, within the ground and 
evapoaration to the atmosphere. It also simulates the movement of solute through soil and 
sediment movement over the land surface. TOPOG is based on the WETZONE program, 
developed by Dr. Emmett O'Loughlin in the mid- 1980. It has undergone several changes 
after that (http://www.per.clw.csiro.au/topog/intro/intro.html). 

TOPOG enables the analysis of three-dimensional landscape through sophisticated digital 
terrain analysis. It is intended for small catchments (up to 10 km^, and generally smaller 
than 1 Icm^). It starts with raw contours and digital elevation model contours. The network 
is divided into cells by series of flow tubes, which diverges or converges, from the ridges 
and valleys on the catchment boundary. Depending upon the shape of the terrain, the 
contours change and hence, the density of the cells can be altered. The standard suggested 
cell size is 20 m by 20 m. The aspect, slope, equinox radiation, the summer radiation, etc. 
are needed to be entered for each cell. The TOPOG calculates the parameters related to 
steady state water balance such as normal wetness index, radiation-weighted index, and 
dynamic water balance parameters such as daily water balance, moisture storage, and water 
table deptli and soil moisture. 
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STREAM 


The Spatial tools for River basins and environment and analysis of management 
(STREAM) is a river basin management instrument. Originally, the RHINEFLOW 
(Kwadijk, 1993; van Deursen & Kwadijk, 1994) was developed to analyse the impact of 
climatic changes on the water resources availability of Rhine basin. STREAM is an 
extension over it and includes the effect of land use changes and coastal dynamics also. It 
is a spatial distributed hydrological model, which uses GIS data to calculate water balance 
in a large river basin. 

STREAM uses land use maps, soil type map, monthly precipitations maps, monthly 
temperature maps and a digital elevation model as inputs. The area under study is divided 
into grids. The water balances for these grids are determined by deriving flow directions 
from DEM and by using Thornthwaite & Mather (1957) equation. During calculation, it is 
able to accounts for human interventions in tenns of land use changes and external 
influences such as climatic changes. The model provides potential and actual 
evapotranspiration, soil water regime, aridity index, runoff regime and snow cover as an 
output. The model has been applied to many large basins, such as Ganges-Brahmaputra- 
Meghna, Rhine, Zambezi and Yangtze basins (Aerts, Kriek and Schepel, 1999). 

The development of the model was done within the framework of the Indian - Netherlands 
water and coastal management cooperative programme. Aerts (2002) describes the 
objective of the work as to develop an instrument which can be used to give insight in the 
long term impacts of land use change, climate change, river basin management, population 
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pressure, economic development and on the future water demand and water availability in 
the Krishna basin of size 258,948 km^. The resolution of the data set is of 1 square 
kilometer. The potential evapotranspiration was calculated using the Priestly - Taylor 
method using monthly global datasets for air temperatures, radiation, albedo, cloudiness, 
and elevation. The model produces water availability and aridity maps. 

Surface Water Model 

Utrecht University developed a surface water model for Orinoco river basin in South 
America, to provide Governments, Non govenunent organisations (NGO's) and other 
stakeholders with comprehensible tools for rough impact assessment of planned human 
interventions in the river basin on surface water runoff (Schot et. al., 2001). The Orinoco 
River is the third largest river in the world having draining area of 830,000 km^. Due to the 
unavailability of large data sets, the model developed by Van Deursen and Kwadijk (1994) 
was used; which requires little input data and can be applied to large river basins like 
Rhine and Ganges-Bralimputra. 

In the model, the water balance consists of three storages, snow, soil water and 
groundwater. These storages are connected by six fluxes, precipitation, snowmelt, and 
evapotranspiration, groundwater recharge, rapid runoff and delayed runoff. There are five 
controls on these fluxes, which are temperature, maximum water storage capacity of the 
soil, crop coefficient, separation coefficient and recession coefficient. The model is based 
on mainly three assumptions. First is that water balance on a location only depends on 
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local conditions from the previous time step and there is no contribution from the other 
locations. Second assumption is that all runoff is assumed to leave the basin within one 
model time step. Third assumption says the precipitation falling at the end of a time step is 
considered as discharge for that time step. 

Precipitation, temperature, maximum water storage capacity of the soil, crop coefficient, 
separation coefficient, and recession coefficient are used as input to calculate the water 
balance per cell of the raster GIS. The flow path of the runoff is decided using digital 
elevation model of the catchment. The runoff volume is calculated by multiplying the cell 
area with the runoff from each cell. By routing the accumulated runoff towards the 
outflow point, local drainage directions were determined. The output from the model is a 
time series of maps showing spatial distribution of flux or storage. The effect of human 
interventions on river discharge was studied by creating hypothetical dam and water 
diversion. 


The WBM Model 

The Water Balance Model (WBM) is a physically based, spatially distributed parameter, 
conceptual catchment hydrological model. It was developed by Prof. Manfred W. 
Ostrowski in 1992 at the Institute for Hydraulics and Water Resources Engineering, 
Darmstadt University of Technology, Germany. It has undergone several modifications 
after that by a group of researchers under his guidance. Since WBM has been employed in 
the present research work, the detailed description of the model is presented in Chapter 3. 



Chapter 3 


The WBM 


3.1 General 


This chapter presents the detailed description of WBM model and its component. The 
WBM uses raster element which represents a column having four layers. The four layers 
are; 


1 . Infiltration layer, 

2. Root layer, 

3. Transport layer, 

4. Ground water layer. 


In each raster element, a detailed simulation of the soil - vegetation - atmosphere - 
transfer - system (SVAT) has been applied which accounts for separation of precipitation 
into storage processes in a vertical column and runoff generation. Each element can be 
partly vegetated as well as partly urbanised. The vertical column with all layers is shown 
in Figure 3.1. 
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Figure 3.1: Vertical column considered for a raster element 

A macro pore system is defined parallel to the first two layers i.e. infiltration layer and root 
layer. The Infiltration layer accounts for infiltration and evaporation. The root layer 
extends up to the maximum rooting depth of the plants and accounts for transpiration. In 
transport layer, water percolates to the ground water layer. The transport processes 
between adjacent elements through these layers are well explained in Figure 3.2. In case of 
high ground water level, capillary rise takes place in the root zone. Backwater effects have 
to be considered for lower two layers. Once the precipitation saturates the macro pore 
system and infiltration layer, surface runoff takes place. The lateral transport of the water is 
always unidirectional and follows the cascade system determined by surface topography. 

For each elementary surface, the water balance is illustrataed with all relevant hydrological 
processes, so that, at each time and at each place, detailed evaluation can be possible. 
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3.2 Coupling with GIS 

Topography based modelling started since 1990s. To model the transport process within 
the catchment, topographic information is necessary. Digital elevation model (DEM) 
pi'ovides such information automatically about the catchment area borders, catchment 
surfaces and direction of flow. 

The Triangulated Irregular Network (TIN) model is a significant alternative to the regular 
raster of a DEM, and has been adopted in numerous GISs and automated mapping and 
contouring packages. In case of TINs or contour representations, sometimes, it is not 
possible to illustrate exactly the salient points like steep embankment of the catchment etp. 
Hence, raster format for the elementary surface is required which provides simple data 
storage and processing. 

To administer the spatial distribution of the land use and parameters, which influences 
different processes, WBM model has been coupled with GIS. IDRISI is chosen as GIS 
program because of its comparatively simple data base management and user friendly 
routine. IDRISI was developed at Clark University. The model describing and model 
steering parameters are made available by providing different ASCII files to the program. 
The data structure and flow chart of the complete system of “IDRISI - WBM - FLIESS” is 
shown in Figure 3.2. The ability of IDRISI towards analysis of a digital elevation mfip 
(DEM) is limited, FLIESS was developed to performe a topopraphic analysis, and produce 
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ecaxtly the output that is needed for WBRM. FLIESS is a topographic analysis tool, that 
finds the flow paths through a DEM by the D8-Method with some additional 
difications, delineates watersheds, gives the slope angles between raster grids and 

produces output for WBRM. 





Figure 3.2: Flow chart of “IDRISI - WBM - FLIESS” 


Note: Flachendaten means flat data, Eingabe means Input. Ausgabe means output. Relationale Datenstruktur 
means relational data structure. 


3.3 Input and output of the WBM 

To caloulate the water balance, the following inputs are required at the desired time 
interval for which runoff is expected fern the model. To mention the units of following 
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parameters, one hour is used as time interval. The WBM is flexible to incorporate any 
time interval. So the units will be changed accordingly. 

i. Precipitation (mm) 

ii. Temperature (°C) 

iii. Evapotranspiration (mm/hr) 

iv. Digital elevation model 

V. Classification of soil in the catchment 

vi. Type of soil with their depth at a particular raster element 

vii. Dimensions of the channel present in the catchment 

viii. Land use characteristics. 

The output from the model contains the values of following at an assigned time interval. , 

i. Through rainfall (mm) 

ii. Interception evaporation and Interception (mm) 

iii. Soil evaporation (mm) 

iv. Transpiration (mm) 

■ v. Evapotranspiration (mm) 

vi. Soil moisture profile for different layers of soil 

vii. Surface runoff (m^/s) 

viii. Base flow and groundwater recharge (m^/s) 

ix. Interflow (mVs) 

X. Macropore flow (mVs) 
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3.4 Model components 

Interception 

When it rains over a catchment, not all the precipitation falls directly onto the ground. 
Before it reaches the ground, a part of it may be caught by the vegetation and subsequently 
evaporated. The volume of water so caught is called interception. Water may be retained 
by the vegetation as surface storage and returned to the atmosphere by evaporation. It can 
drip off the plant leaves to join the ground surface or the surface flow. The rainwater may 
run along the leaves and branches and down the stem to reach the ground surface. Change 
in Interception storage can be described using following function: . 

at 

Where, p: Passage coefficient; N: Precipitation; D; Drainage; E: Interception evaporation; 
S: Storage 

Discharge component 

Topographic information derived from DEM provides the basis for simulation of lateral 
water transportation. It has mainly three components - (1) the surface discharge, (2) the 
macropore storage and (3) the interflow. Each raster element receives water from its 
neighbours and passes on to its successors as shown in Figure 3.3. The arrows point to the 


successor. 
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Figure 3.3: Cascade model of raster elements. 


Surface discharge 

Land surface discharge moves as layer discharge in small rivulets until a larger stream 
forms. Each raster element is covered with complete water film which flows off over the 
entire width of an edge as shown in Figure 3.4. The discharge behavior is described by 
continuity equation and Manning Strickler formula. 
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Figure 3.4: Schematic representation of a raster element 


dt 



a 

A 


Where Nefr = Effective precipitation; Qa = Surface discharge; Qz= Lateral discharge; 
A = SurLice area of the raster element; h = height of water film 


Where Vm = Velocity of the flow; Ks = Roughness factor; i = download gradient as shown 
in Figure 3. 4. 


Qa=V„hB 


Where B is the width of the discharge. Total discharge is counted on all elements of the 
catchment area whose surface exceeds a certain threshold value. 
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Macropore discharge 

Capillary forces play a crucial role in macropore discharge. The description of the 
macropore system is very difficult to sieze. The size of the macropore system is expressed 
as percentage share of the total pore volume (TPV). Rain water, which cannot infiltrate 
into the soil matrix, first fills the macropore volume and when the capacity of the 
macropore is completely filled, it adds to the surface discharge. 

Soil Moisture Simulation 

A soil column consists of three homogeneous layers. Transportation of water from these 
layers is shown in Figure 3.5. 



Transpiration 


Infiltration 

1 

A 
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k. 
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rise 

1 

Groundwater 

recharge 

Ground water layer 

— 


Interflow 
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Figure 3.5: Transportation of rain water through Soil columns 

The computation of water balance at time t is calculated from the following relation: 

~ = I{t) + K{t) - Pit) - Ifit) - E{t) - Tit) 
at 

Where M = Soil Moisture (mm); 

I(t) = rate of actual Infiltration (mm/h); 

K(t) = rate of Capillary rise (mm/h); 

P(t) = rate of actual exfiltration (mm/h); 

If(t) = rate of actual Interflow (mm/h); 
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E(t) = rate of actual Evaporation (mm/h); 

T(t) = rate of actual Transpiration (mm/h). 

While calculating the soil moisture balance, it is to be noted that in layer 1 , capillary term 
will be zero; in layer 2, Evaporation term will be zero; in layer 3, capillary, interflow, 
evaporation and transpiration term will be zero. 

Infiltration 


When water is applied to the surface of a soil, a part of it seeps into the soil. This 
movement of water through the soil surface is known as infiltration. The rate of infiltration 
has been calculated using following relationship: 


UFc-ejt) 


fi.FC-fi.FC 


+ / 


end 


Where, 0 = actual soil moisture [mm]; Iraax= rate of maximum Infiltration [mm/h]; 

lend = final rate of infiltration [mm/h]; f, = starting empirical factor; f, = final empirical 

factor; PC is field capacity. 


For 0 < f,*FC, I(t) = Imax; For 0 > £*FC, I (t) = lend. 
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Percolation 


To calculate the percolation, following equation has been used: 


P(t)^Kf. 


e{f)-f^.Fc 

UFC-f,.FC 


Where Kf — Saturated hydraulic conductivity [mm/h]; a = Soil type specific Exponent 


Interflow 

A part of the precipitation that infilters moves laterally through upper crusts of the soil and 
returns to the surface at some location away from the point of entry into the soil. This 
component of runoff is known as interflow. 

To calculate the interflow, following equation has been used: 


7(0 = sin 


eit)-f.FC 

YUFC-f.FC^ 


Where P = angle of inclination. 


For0<f,*FC, P(t)/I(t) = 0; 

For 0 > fi*FC, P(t)/I(t) = sinp * Kf . 
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Capillary rise 


The soil water is broadly classified into two categories: (1) Free water, and (2) Held water. 
Free water moves in the pores of the soil under the influence of gravity. The held water is 
retained in the pores of the soil, and it cannot move under the influence of gravitatioriial 
force. Maximum rise of capillary has been calculated using following equation: 


m=i- 


0(t)-f,.FC 

f,.FC-f,.FC 



Where, Kmax = Maximum rate of capillary rise [mm/h] 


For 0 < f,*FC, I(t) = Kmax; 
For 0 > fi*FC, I(t) = 0. 


Evapotranspiration 


Evapotranspiration is the sum of evaporation and Transpiration. Evaporation comes firom 
the highest layer (Infiltration layer). The maximum evaporation is reached with field 
capacity. Transpiration takes place from second layer and depends upon the available soil 
moisture. The following relations have been used to calculate these variables. 




m 

FC 
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T{t) = [Ept-Em 


e{t)-WP 

FC-WP 


With the knowledge of model components and required input data, the next chapters deals 
with the data processing to generate suitable data to apply the model to a catchment. 



Chapter 4 


Study area and data processing 


4.1 General 

In the present study, WBM model has been used to simulate various hydrological 
processes in Tarrawarra catchment, Australia. Due to sufficient online available data, the 
Tarrawarra catchment was selected which is located near Melbourne in Southern Victoria, 
Australia. In this chapter, first the description of the catchment along with available data 
has been presented; which will be followed by a brief description of the required inputs 
and obtained outputs from the WBM model. 

4.2 Description of the catchment 

The area of Tarrawarra catchment (37°39' south, 145°26' east) is 10.5 hectares. The 
catchment has smoothly undulating terrain with no perennial streams and no channels. It 
has few 10 m cypress trees, 5 m mixed Australian native trees and two Eucalyptus trees. 
The vegetation cover is greater than 90%. The land is used for cattle grazing. A road 
follows along the northern boundary. The runoff from this road does not contribute to the 
catchment. The aerial photograph of the catchment looking south east is shown in Figure 


4.1. 
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The climate is temperate with mean annual rainfall of 820 mm and potential 
evapotranspiration of 830 mm. The minimum and maximum mean temperatures in 
summer are 9°C and 24°C, respectively; and those for winter are 3°C and lO^C, 
respectively. I 

The soil in the catchment consists of 3 layers. At the upper slope of the catchment, it has 
typically 20 cm thick silty loam or silty clay loam A1 horizon, a silty clay loam A2 horizon 
and heavy clay B horizon. The midslope soils have a loam clay loam A horizon of 20 cm 
thick, which gradually transitions into a silty clay loam B horizon. The surface soil cracks 
during dry periods. The definition of different horizons can be found in Appendix B. 



Figure 4.1: Aerial photo of Tarrawarra catchment looking south east 
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4.3 The Tarrawarra data set 

The Tarrawarra data set are available on Dr. Andrew Western's home page, who is Senior 
Research Fellow in the Cooperative Research Center for Catchment Hydrology and Center 
for Environmental Applied Hydrology at the Department of Civil and Environmental 
Engineering at the University of Melbourne, Australia. It can be directly dovraloadable 
from the following link: 

http://www.civag.unimelb.edu.au/~westem/tarrawarra/tarrawarra.html . 

The available data has been classified into 8 categories. The content of each file is shortly 
described below. 

Topographic Data 

A detailed topographic survey on a 10 m by 10 m paced grid was conducted within the 
local coordinate system using Total Stations. Two grided digital elevation models (DEM) 
of 5 m grid have been developed with Tarrawarra coordinate system and Universal 
Transverse Mercator (Zone 55). 

Surface Soil Moisture Data 

Soil Moisture maps of 10 m by 20 m grid have been collected using time domain 
reflectometry (TDR) over the top 30 cm soil layer. To ensure the accuracy of measured 
data, TDR data were calibrated with gravimetric measurement of soil moisture under some 
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typical field conditions. The data file contains date, time, point coordinates, dielectric 
constant and moisture data. 

Soil Moisture Profile 

Soil moisture profiles for 20 sites were measured using Campbell Pacific Nuclear 503 DR 
Neutron moisture meter over a period of 2 years. Measurements were taken by access tube 
at 15, 30, 45, 60, 90, 120 and 150 cm, depending upon the depth to bedrock. The site co- 
ordinates, depth to bedrock and soil profile description are available in 20 files. 

Piezometer Data 

Seventy four piezometers were installed to depths of between 0.5 m and 1.1 m to measure 
shallow water table elevations. Date, time and depth to water table from natural surface 
for 74 points are listed in a file. 

Soil data 

Saturated hydraulic conductivity, particle size data and soil profile data are available for 
the catchment. The soil profile data consists of point coordinates, the depth to the bottom 
of the A horizon (cm), the depth to the bottom of the B1 horizon (cm), the texture category 
of the B1 horizon, the depth to the bottom of the B2 horizon (cm), the texture category of 


the B2 horizon. 
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Metereological and Runoff Data 


Meterological data from an automatic climate station and 5 collecting rain gauges are 
available at 6 min. interval for 3 consecutive years (1995 to 1997). The Meterological data 
consists of dry bulb air temperature, wet bulb air temperature, rainfall, average wind speed, 
maximum wind speed, wind direction, global radiation and net radiation. A separate file 
with daily summary of minimum, maximum and mean values of temperature, daily 
accumulation of rainfall and radiation variables, wind - run and maximum wind speed are 
available. 
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Figure 4.2: Tarrawarra catchment with contours lines and meteorological stations 

To measure the surface mnoff, a flume has been installed near the outlet of the catchment. 
The positions of climate stations, rain gauges and flume are shown in Figure 4.2. Water 
levels in the flume are measured every minute and are saved if the change in water level 
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exceeds 1.5 mm from the previous measurement. This means that if the change is less than 
1.5 lum, then it is not recorded and it is assumed to be constant. With this logic, when the 
data are not available (means not recorded) then the actual value of flow is within 1.5 mm 
of the last recorded flow value. Therefore, the missing flow values can be assumed 
constant equal to the last recorded value. 

jl 

Surface Roughness and Vegetation data 

Vegetation was sampled from Im by Im plot. The vegetation consists predominantly of 
perennial and annual introduced grass and clover pasture species typical of Australian 
legion. Sampling coordinates, water content and biomass in g/m^ have been provided in a 
file. Surface roughness was measured using pin profiler having 41 pins at a spacing of 25 
mm. Six profiles were measured at each site: two across the slope; two up the slope; and 
two at 45° to the slope. 

Photographs and Maps of the catchment 

Contour map, moisture maps, aerial photographs looking catchment from different 
directions, photographs of different instruments used in data acquisition was available. 
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4.4 Input files used in the model 


The input files with their short description are presented below. 


BNR_ATAR 


The file contains information regarding region showing similar soil type in the catchment. 
The entire catchment has been divided into 20 different regions by assuming 20 different 
points and applying Theisen polygon method. 


DHM_ATAR 



Figure 4.3: Digital elevation model of Tarrawarra catchment 
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GCA_ATAR 

This contains information about the number of grid cells that are contributing flow to the 
actual grid cell. 


GEF_ATAR 


This provides information about the slope of a grid cell to its consecutive grid cell as an 
input to the model. 

LNR_ATAR 

This file contains infonnation about the distribution of land uses. 


NAC_ATAR 


This conveys about the flow pattern to the model. This file contains information that how 
many grid cells are contained in each grid cell. 

PZR_ATAR 


This file contains information regarding the regions for which soil moisture profile is 
desired after simulation. 
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Weiras.AAW 

The duration of simulation run with time interval is defined in this file. 

Weiras.BOA 

The type of soil within the catchment with their welting point, field capacity, saturated 
hydraulic conductivity and maximum infiltration is assigned in this file. 

Weiras.BOD 

The depth of soil layers with their type for all the categorised regions within the catchment 
is assigned here. 

Weiras.EFE 

In this file, user can define the information about number of different land uses with 
number of regions with similar soil types in the catchment. For example, in the present 
study, the entire catchment is considered having single type of land usage and the soil 
within the catchment has been divided into 20 regions. 
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Weiras.EFT 

The threshold values necessary to form a channel is assigned here. 

Weiras.EZG ; 

This file associates the input time series of a catchment. 

Weiras.FIL 

The name of input rain, temperature and evapotranspiration file with their unit and time 
interval in which these data are available for the model. 

WeirasJGG 

This file contains the factors by which the values of parameters assigned in Weiras.LNZ 
will be multiplied and results into the new values of parameters for different month. 

WeirasXNZ 

The root depth, vegetation cover and its seasonal variations, leaf area index and its 
seasonal variation, percentage of unproductive land, surface roughness value and Haude 
monthly coefficient for evapotranspiration is assigned in this file. 
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Weiras.TRS 

The geometry of channel is specified here in this file. 

4.5 Generation of input files from available data 

The required input files along with its preparation from available data has been discussed 
in this section. 

4.5.1 Soil classification 


The soil profile data of catchment was available. The available data files looks like as 
follows. 




A horizon 

B1 horizon 

B1 horizon 

B2 horizon 

B2 horizon 

X 

y 

depth to 

depth to 

texture 

depth to 

texture 

coord 

coord 

bottom 

bottom 

category 

bottom 

category 

m 

m 

cm 



cm 


805 

950 

24 

39 

clay 

>72 

clay 

805 

930 

24 

46 

clay 

>75 

clay 

805 

910 

16 

51 

clay 

>74 

clay 


B 1 horizon texture was taken as reference to categorise the available data into three types 
of soil - clay, silt and silt clay. The depth of layers was calculated from here. Twenty 
points were chosen out of this data set and used to create BNR_ATAR file. The generated 
output BNR file is shown in Figure 4.4. 
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Figure 4.4: Regions with different soil types 

The soil classification with information about different layers has been provided in 
Weiras.BOD file. The input Weiras.BOD for the present study is shown below. 


Boden C'^,BOD) 


ID 

mak 

anzsch dl 

boal 

Boden 
d2 boa2 

-Kenngroessen 

d3 boa3 d4 boa4 d5 boaS d6 boa6 

GW 

Stand 

Beschreibung 


mm 


m 

- 

m 

- 


m 


- 






1 

A 

1 

1 

Y 
I 

A 

1 

1 

Y 
1 

A 

1 

1 

Y 
1 

A 

1 

1 

Y 
1 

A 

1 

i 

Y 

A 

1 

» 

Y 
1 

A 

1 

1 

Y 
1 

A 

1 

1 

Y 
1 

A 

1 

1 

Y 

< — > 



— -> 

1 


3 

0.21 

2 

0.18 

3 

0.36 3 


Boden 

I 


2 


3 

0.19 

2 

0.24 

3 

0.32 3 


Boden 

II 


3 


3 

0.26 

2 

0.27 

2 

0.22 3 


Boden 

III 


4 


3 

0.25 

2 

0.22 

2 

0.28 3 


Boden 

IV 


5 


3 

0.19 

2 

0,24 

2 

0.32 3 


Boden 

V 


6 


3 

0.19 

2 

0.36 

2 

0.20 1 


Boden 

VI 


7 


3 

0.19 

2 

0.29 

2 

0.27 1 


Boden 

VII 


8 


3 

0.25 

2 

0.30 

2 

0.20 1 


Boden 

VIII 


9 


3 

0.22 

2 

0.27 

2 

0.27 1 


Boden 

IX 


10 


3 

0.24 

2 

0.14 

2 

0.37 2 


Boden 

X 


11 


3 

0.21 

2 

0.22 

1 

0.32 3 


Boden 

XI 


12 


3 

0.22 

2 

0.22 

1 

0.31 3 


Boden 

XII 


13 


3 

0.17 

2 

0.15 

1 

0.43 3 


Boden 

XIII 


14 


2 

0.25 

2 

0.50 

3 



Boden 

XIV 


15 


2 

0.22 

2 

0.53 

3 



Boden 

XV 


16 


2 

0.26 

2 

0.49 

2 



Boden 

XVI 


17 


2 

0.24 

2 

0.51 

2 



Boden 

XVII 


18 


2 

0.30 

2 

0.45 

2 



Boden 

XVIII 


19 


2 

0.29 

2 

0.46 

2 



Boden 

XIX 


. 


2 

0.23 

2 

0.52 

1 



Boden 

XX 



Note: Boden means soil; Beschreibung means description; boal, boa2 indicates the types of soil and dl, d2 
are depths of different soil layers. 
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4.5.2 Rainfall 

There were five rain gauges installed in the catchment to measure rainfall. The rainfall 
data along with other meteorological data at 6 min interval was available for three years 

I 

1995, 1996 and 1997. The daily summary of the meteorological data was also available in 

a separate file. To verify the uniformity of rainfall throughout the catchment, it was 

necessary to compare the available rainfall data for all five rain gauges. The time series 

graph showing rainfall from different rain gauges are shown in Figure 4.5. This figure 

clearly indicates that the rainfall can be assumed to be uniform for whole catchment. 

$ 

Rainfall values in the required format were prepared to provide as an input to the model. 



Figure 4.5: Time series graph of available rainfall 



42 


4.5.3 Temperature 


The diy bulb air temperature, wet bulb air temperature, surface air temperature were 
available. The surface air temperature has been used in the model. 

4.5.4 Evapotranspiration 

The meteorological data were used to calculate evapotranspiration using Permmaim 
Monteith Equation which has been described later in this chapter. 


4.5.5 Runoff 


Instantaneous surface runoff values for three years were available. The available data set of 


observed surface runoff looks as follows. 


Units: runoff is in mm/h 
20 / 10/95 8 : 30:00 0.0194 
20 / 10/95 8 : 49:00 0.0172 
20 / 10/95 9 : 03:00 0.0148 
20 / 10/95 9 : 09:00 0.0127 
20 / 10/95 9 : 35:00 0.0105 
20 / 10/95 10 : 00:00 0.0081 
20 / 10/95 10 : 43:00 0.0069 
20 / 10/95 11 : 28:00 0.0059 


In the present case the defined time interval was 6 min. So, the obtained surface runoff 


from model was at 6 min interval only. To compare the time series graph, it was necessary 


to extract data of observed surface runoff at 6 min. interval. From the available data set, the 


values which were lying at time step which are multiple of 6 min. were only captured and 
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others were discarded, e.g. out of data set shown above, the captured surface runoff values 
will be at 8:30 and 10:00, because only these times are falling at a multiple of 6 minute. 
The obsei'ved surface runoff directly drawn from available data and the extracted observed 
surface runoff are shown in Figure 4.6a and Figure 4.6b, respectively. By observing these 
two figures, it is clear that the information has not been lost due to this methodology to 
extract the available data as per the required format. 



26.3.96 15.5.96 4.7.96 23.8.96 12.10.96 

Date 

Figure 4.6a: Time series from the available observed data 




Runoff (Cumecs) 
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20000 30000 40000 50000 60000 70000 

Time from 1/1/96 to 31/12/96 

Figure 4.6b: Time series from extracted data 


4.6 Output files obtained from the model 

There are six folders in which output values are stored after the model run. 

i. OUTWEL 

ii. BODFUN 

iii. OUTSMU 

iv. OUTWBM 
V. OUTPZR 


vi. OUTRST 
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In OUTWEL folder, all the meteorological information is stored. Most of the times, the 
output files stored in this folder is only required for analysis of simulation results. For 
further and deep investigation of results in terms of movement of water in soil layers and 
to the atmosphere, one needs to check files stored in the other folders. 

BODFUN contains file which shows functions generated for infiltration, exfiltration rates 
and percolation with available soil moisture (refer to Figure 5.1 1). The data files stored in 
OUTPZR folder presents values for soil moisture in volume % for all the defined points in 
pzr_atar file and at assigned time interval. OUTRST contains soil moisture profile in raster 
format as per the information provided in the weiras.EFT file. 

4.7 Data Processing 

The Meteorological data are reported to be available for three years; but actually these are 
available for the following period with lot of missing values. 

Year 1995: 10/8/95 to 31/12/95 
Year 1996; 1/1/96 to 31/12/96 
Year 1997; 1/1/97 to 17/11/97 

By deep investigation of the available data, the following five periods were determined, in 
which the continuous 6 minutely meteorological data are available. The periods are 


mentioned below: 
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Period 1 : 12:06 on 26/6/96 to 23:54 on 3/7/96; 

Period 2: 00:00 on 17/7/96 to 21:24 on 1 1/8/96; 

Period 3: 00:00 on 24/8/96 to 12:06 on 31/8/96; 

Period 4: 09:30 on 2/9/96 to 04:06 on 16/9/96; 

Period 5: 00:00 on 30/9/96 to 12:48 on 8/10/96. 

The runoff values in this period were almost continuous except some periods. The missing 
runoff data were filled using linear interpolation. To run the model, evapotranspiration 
values at 6 minute interval for these periods were required. The widely used Penman - 
Monteith equation calculates the evapotranspiration values only at the daily basis. To 
apply the Penman - Monteith equation on an hourly or shorter time scale, the equation 
needs to be adjusted. Smith (1992) suggests the Penman-Monteith equation for hourly 
time step as: 


0.408A(i?„-Cr) + y. 


37.5 


ET,=- 


T,,+213 


U2(.e\TJ-e,) 


A + y(l + 0.34 m2 ) 


ETo is the reference evapotranspiration (mm hour"’); 

Rn is the net radiation at the grass surface (MJ m'^ hour'*); 

G is the soil heat flux density (MJ m'^ hour'); 

Thr is the mean hourly air temperature (°C); 

A is the saturation slope vapour pressure curve at Thr (kPa °C"'); 
y is the psychrometric constant (kPa °C''); 



47 

e°(Thr) is the saturation vapour pressure at air temperature Thr (kPa); 

Ca is the average hourly actual vapour pressure (kPa); 

U 2 is the average hourly wind speed (m s‘‘). 

Evapotranspiration values at 6 minute interval can be calculated using above equation by 

/ 

changing factor from 37.5 to 3.75'in the second term of numerator. The other variables also 
need to be provided at 6 minute interval in order to get result at 6 minute interval. 

So, the equation becomes, 


0.408A(i?„ -G) + y. 


3.75 


^ro=- 


T,,+273 




A + y(l + 0 . 34 ^ 2 ) 


Given relative humidity measurements, the actual vapour pressure is determined as: 


e 


a 




RH,r 

100 


Where, eg is the average actual vapour pressure (kPa); e°(Thr) is the saturation vapour 
pressure at air temperature Thr (kPa); RHhr is the average hourly relative humidity (%). 


(r) = 0.6108.exp 


17.27 
r + 237.3 


Where, e°(T) is the saturation vapour pressure at air temperature T (kPa) and T is the air 
temperature (°C). 
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4098- 


A = - 


0.6 108 exp 


(17.27*r) ' 
7’ + 237.3 


(r + 237.3)^ 


Where A is saturation slope vapour pressure curve at T (kPa 


G = 0.1*Rn during daylight periods; G = 0.5*Rn during night time periods. 


4.8 Calibration and Validation data set 

The minimum, maximum and standard deviation of relative humidity, ground temperature, 
rainfell, wind speed, runoff and net radiation for above mentioned 5 different periods are 
tabulated below. 



RH 

GRND 

RAIN 

WIND 

Runoff 

NetRad. 


% 

TEMP 

6inin 

6min 

M3/sec 



* 

Mean 

total 

mean 


mean 


— 

calc 



inm 



W/m^ 

1 Period 1: 12:06 on 26/6/96 to 23:54 on 3/7/96 (Total data: 1799) I 

Max 

97 

15.3 

1.2 

23.8 

WBEISaM 

337.8 

Min 

87 

2.3 

0 

0.4 


-80.2 

Average 

91.93 

9.29 

0.01 

6.12 

0.0013 

11.90 

Std 

deviation 

1,69 

2.45 

0.07 

4.65 

().()026 

71.46 

Period 2: 00:00 on 17/6/96 to 21:24 on ll/8/96(Total data: 6211) 

Max 

94 

16.5 

1 

38.3 

0.0511 

427.8 

Min 

38.00 

1.50 



0.0000 

-81.50 

Average 

80.37 

9.18 

0.01 

7.85 

().0()18 

14.96 

Std dev 

8.08 

2.66 

0.06 

5.51 


79.67 

Period 3: 00:00 on 24/8/96 to 12:06 on 31/8/96(Total data: 1802) 

Max 

93 

19.1 

0.6 

20.8 

0.0104 

471.7 

Min 

52 

2.9 

0 

0.4 

0.0000 

-79.3 

Average 

77.67 

10.53 

0.01 

J^oo 

0.0012 

37.86 

St dev 

10.62 

3.34 

0.06 

3.60 

0.0019 

109.35 














































49 



Period 4: 09:30 on 2/9/96 to 4:06 on 16/9/96(Total data: 3306) 

Max 

91 

19,1 

2.2 

41.6 

0.0541 

578.6 

Min 

48.00 

3.40 

0.00 

0.40 

0.0000 

-75.70 

Average 

74.66 

9.98 

0.02 

10.31 

0.0027 

43.55 

St dev 

9.11 

2,76 

0.09 

6.49 

0.0048 

115.34 

Period 5: 00:00 on 30/9/96 to 12:48 on 8/10/96( Total d: 

ata: 2047) 


Max 

92 

23 

1.4 

31.1 

0.053 1 

747.6 

Min 

42.00 

5.60 

0.00 

0.40 

0.0000 

-71.40 

Average 

73.25 

12.22' 

0.02 

8.41 

0.0017 

78.07 . 

St dev 

11.36 

3.39 

0.08 

5.53 

0.0046 

162.32 

Overall (Total data: 15165) 

Max 

97 

23 

2,2 

41.6 

0.0543 

747.6 

Min 

38 

1,5 

0 

0.4 

0.0000 

-81.5 

Average 

79.22 

9,94 

O.Ol- 

7.93 

0.006 1 

32.07 

Std Dev 

10.26 

3.03 

0.07 

5.71 

0.0102 

107.55 


By investigating the above statistics, two periods were selected for further analysis, (1) 
Period“4 (total data points 3306) for calibration and (2) Period-2 (total data points 6211) 
for validation. 


It seems that Period 1 and period 3 are of limited range only, so these periods were not 
considered for further study. The time series plot for calibration and validation periods are 
shown in Figure 4.7a and Figure 4.7b respectively. 














Runoff (Cumecs) 



Figure 4.7a: Time series plot for calibration period 



Figure 4.7b: Time series plot for validation period 
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The calibration period-4 was further classified into 3 sets on the basis of low peak flow 
(denoted by period 4L), medium peak flow (denoted by period 4M) and high peak flow 
(denoted by period 4H). The periods are as follows: 

Period 4L: 9:30 on 2/9/96 (1) to 1 1:24 on 4/9/96 (500); / 

Period 4M: 1 8:24 on 1 1/9/96 (2250) to 10:24 on 13/9/96 (2650); 

Period 4H: 15:24 on 13/9/96 (2700) to 21:30 on 14/9/96 (3001). 

The numbers in the bracket corresponds to the number used to show time at 6 min Interval 
on the x-axis of the graphs shown in Figure 4.7a and Figure 4.7b. 

Similarly the validation period were also divided into 3 categories. The periods are as 

i 

follows: 

Period 2L: 10:12 on 24/7/96 (1783) to 13:24 on 27/7/96 (2535) 

Period 2M: 16:54 on 20/7/96 (890) to 17:42 on 21/7/96 (1 138) 

Period 2H: 4:24 on 29/7/96 (2925) to 20:24 on 30/7/96 (3325) 





Chapter 5 


Parameter estimation 


5.1 General 

Increasingly sophisticated computer-based hydrologic models are currently used in a 
variety of tasks. These hydrologic models must be calibrated to be useful for the solution 
of practical problems. Through calibration, we estimate values for the model parameters, 
which enable the model to closely match the behaviour of represerited real system 
(watershed). It is generally believed that appropriate values, for certain model parameters, 
can be determined through direct measurements conducted on the real system. However, in 
many situations, the model parameters are conceptual representations of abstract physical 
characteristics and must be adjusted through a trial and error process that matches the 
model response to historical input output data. 

There is no universal systematic procedure for model calibration. In this section, a 
particular sequence of calibration procedure is discussed, which has been adopted in the 
present work. The idea behind calibration of the model is to estimate the following seven 
parameters: 



i. Infiltration rate (Kf) of silt clay in (mm/hr) 

ii. Kf of silt in (mm/lir) 

iii. Kf of clay in (mm/hr) 

iv. Maximum infiltration rate of silt clay in ( mm/ hr) 

V. Maximum infiltration rate of silt in (mm/hr) 

I! 

vi. Maximum infiltration rate of clay in (mm/hr) 

vii. Roughness (Kst) 

The values of these parameters should be chosen in such a fashion that the output runoff 
from the model should approximately match with the observed runoff The optimum 
values of these parameters were desired. Root mean square error (RMSE), normalised root 
mean square error (NRMSE), Nash Sutcliffe Efficiency (E), Coefficient of efficiency (El), 
Correlation coefficient (R), average absolute relative error (AARE) and sum squared error 
(SSE) were also calculated while calibration and validation of the model. These error 
statistics can be calculated using following equations. 


RMSE = 



'y'j (Sobs Qcal ) 


( 1 ) 


NRMSE = 



( 2 ) 
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El- 


z 

Qobs Qcal 

z 

Qobs Qobs 


(4) 


R = 


t 

^ 1 Qobs 

>(j 

2ca/ ~ Qcal )). 

|2 

I 

j i^obs Qobs ) 

♦I 

J i^cal Qcal 

J 


(5) 


AARE = —Y 


Qcoi-Q, 


obs 


Q. 


obs 


‘100 


( 6 ) 


sse = Y^{q„,-q,J 


MBE = 100* 

/ j Qobs 


0 ) 

( 8 ) 


The step-wise procedure adopted in calibration of the model is described below. 

First of all, select data set for which maximum numbers of continuous runoff, rainfall, 
temperature and meteorological data to calculate evapotranspiration are available at 6 
minute interval. Then divide the data set for smaller periods so that model run should not 
take much time. To begin with, assign a guess value for the above mentioned 7 parameters 
and run the model. Draw time series graph with observed and calculated data and calculate 
error statistics. Check the sensitivity of all the parameters and accordingly decide the 
value of the parameters. Apply the selected combination of parameters for another period 
and check the time series and error statistics. Once the best match in time series is 
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obtained the combination of parameters will be declared as optimal from trial and error 
procedure. 

5,2 Calibration using 1995 data set 

First the model calibration was tried with the available data for 1995. There were two 
different data sets in which the data were continuous. The results obtained while 
calibrating the model using these data set are presented below. 

5.2.1 Trial 1 

The period chosen for this case study was from 8:30 on 20/10/95 to 7:00 on 22/10/95. 
Arbitraiy combinations of seven parameters were tried. For each combination, time series 
graph for observed and predicted runoff were drawn in one sheet. Fifteen different 
combinations of values were tried. The best match in time series graph was observed for 
the values mentioned in Table 5.1. The time series plot for this combination is shown in 
Figure 5.1. 


Table 5.1: Best combinations achieved from arbitrary selection 


Kf (mm/hr) 

4.5 

4.5 

0.05 


Maximum Infiltration 
(mm/hr) 

4.6 

4.6 

0.06 




3.0 
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Figure 5.1: Time series graph for trial 1 

The predicted runoff patterns were similar to the observed runoff. For better observation, it 
was necessary to increase the duration from 2 days to longer periods. Next the sensitivity 
of Ksi was decided to observe. 

The values of Kst were varied from 1 to 9.0 (1.0, 2.0, 2.5, 2.7, 2.8, 3.0, 3.5, 7.0, 9.0) to 
achieve match in peak runoff value. Few typical effects are presented in Figure 5.2a to 
Figure 5.2d. For smaller values of Kst, the predicted runoff was not following the patterns. 
For higher values (5.0, 7.0, 9.0), the predicted peak runoff were very high compared to 
observed runoff peak value and after certain point it started decreasing. The similarity in 
time series was observed for Kst values ranges from 2.0 to 3.5. Root mean square error 
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With the Kst value fixed as 3, the peak was almost matching but at initial period, the 
pattern was shifted downward (clear from Figure 5.2b). So, to get the best fit in time series 
plot, the variation in values of other variables was tried by keeping Kst fixed as 3.0. By 
observing the RMSE values for each combination, it was clear that the combination 
selected earlier was the best. 


From this analysis, the following conclusions could be drawn: 

1 ) Best values of parameters are as given in Table 5. 1 

2) The study period should be increased or changed to observe the complete effect and 
to ensure the validity of chosen combination. 










Runoff (cumecs) 
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Figure 5.3e: Time series plot for Kst 9.0 and trial 2 

t 

It was decided to use the data set for which data of observed runoff for longer duration is 
available. From the available data set, the observed runoff data for 1996 was found to be 
suitable. So, the next calibrations were done using runoff values measured for year 1996. 

5.3 Calibration using 1996 data set 


The available data for 1996 were used to calibrate the model. There were three different 
periods in which most of the data were continuous. The results obtained while calibrating 
the model using these data set are presented below. 
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5.3.1 Data of 1996 Set 1 

For this trial, the study period was from 00:00 on 1/5/96 to 12:18 on 18/6/96. In this 
peiiod, only two peaks were there in the plot of time series graphs from observed data set. 
The objective to be attained during calibration was same as mentioned earlier. The earlier 
chosen values of parameters were used to note the effect. The time series graph is shown 
in Figure 5.4. 


0.02 

0.015 

1 ^ 0.01 

o 
a 

3 

0.005 

0 

0 1000 2000 3000 4000 

Time (at 6 min. interval) 

Figure 5.4: Time series plot for 96setl 

Some remarkable points from this plot are time to start of rising limb was matching; time 
at the end of falling limb is not matching; the measured peak runoff values are slightly 
higher than observed value; the base time period of hydrograph obtained from model is 
less than that for the observed hydrograph. 


Obs Fred 

;; 


. i 
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It is clear from Figure 5.4 that the volume of runoff obtained from the observed 
hydrograph is higher than that obtained from the model. So decrease in the infiltration rate 
could be one of the solutions. But this will certainly result into an increase in peak value. 
So, anothei possibility could be to make the ground surface rough, so that it will allow the 
sui face lunolt to stay long and provide sufficient time to percolate. Different combinations 
ol values were tried by decreasing the infiltration rates as well as by increasing the 
roughness values; but the desired similarity in hydrographs could not be achieved. 

Up until now, the infiltration values for clay were kept very small. So, it was necessary to 
check its sensitivity. By increasing the infiltration rate for clay and decreasing the same 
for silt clay and silt form the earlier selected data, it was observed that the peak 
disappeared. So, it could be deduced from here that by decreasing the infiltration rate for 
clay, one can get more number of peaks. The sensitivity of other infiltration values was 
also required to be checked. After several trials of different values of seven parameters, it 
was clear that only by decreasing the infiltration values for silt, the peak has increased. It 
means by increasing these values, the peaks will be decreasing. 

By keeping these observations in mind, 15 different trials were carried to attain the exact 
match in time series graph. The differences in base time period suggests the requirement 
to add base flow to surface runoff; which makes the base time period more. So, it could be 
derived that the runoff from the model does not account for base flow. By observing the 
topography of the catchment, it was very clear that there should be two streams. So 
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piestnct of many peaks as an output from the model could be due to time lag in both the 
surface luuoff. Since the runoff in the field was measured only at the outlet, this reason 
seems to be logical. 


By now, the following necessary requirements were: 

i. There should be some guideline to choose the values for infiltration rate; 

ii. The base flow and topography of the catchment should also be accounted for the 
further analysis; 

iii. By incorporating above two points, a decision space should be fixed, in which the 
optimum values of parameters can be searched; 

iv. While searching the parameters in decision space, there should be some objective 
function to attain. 

5.3,2 Data of 96 Set 2 

The infiltration rates for clay, silt and silty clay were investigated firom the literature. In 
most of the cases, it was observed that the Kf and the Maximum infiltration values for silty 
clay and silt are approximately same, and it is half of that corresponding to clay. The 
maximum infiltration values were chosen 10 % more than the infiltration rate. So, for 
further studies, following relationship was finalised. 
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Type of soil 

Kf (mm/hr) 

Maximum Infiltration 



(mm/hr) 

Silt clay 

al 

l.lal 

Silt 

al 

l.lal 

Clay 

0.5al 

0.55al 


With the help of the above relationships, one needs to choose only one value of al and 
accordingly other values would change. Twenty different trials were carried out to fix the 
K,.jt values. The chosen Kst values varied from 0,5 to 10.0. It was observed that there is no 
change in time series graph for the Kst values from 0.5 to 4.0; the significant changes were 
recorded only for the K,t values 4.0 to 10.0. The RMSE was minimum for the Kst value 
5.0. It was decided to fix Kst as 5.0 for the further analysis. 

For each combination RMSE and NSE were calculated. Fifty different trials were carried 
out for the data set having period 12:00 on 26/8/96 to 00:00 on 17/9/96. The infiltration 
values started from 0.15 mm/hr and increased to 7.50 mm/hr at an interval of 0.15. The 
minimum RMSE and NSE were observed for the combination of values as shown below. 
The time series graph for these values of parameters is as shown in Figure 5.5. 


Kf (mm/hr) 

Maximum Infiltration 

Kst 


(mm/hr) 


5.10 

5.61 


5.10 

5.61 

5.0 

2,55 

2.81 
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11000 13000 15000 17000 

Time (at 6 min. interval) 


Figure 5.5; Time series plot corresponding to minimum error statistics 

By observing the time series graph for all the trials, the best match in observed and 
predicted runoff was observed for the following combination. 

Kf (mm/hr) Maximum Infiltration Kst 

(mm/hr) 

1.05 1.15 

1.05 1.15 5.0 

0.53 0.58 

The tinrc scries graph for this combination is as shown in Figure 5.6. In order to find out 
the reason for the higher peak values, the rainfall was also plotted on the same sheet. The 
peaks obtained from model were very logical. Whenever there is maximum rain or 
consecutive ruins, the peaks were present. Similar was the case with observed data. The 
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peaks were too high and number of peaks was nearly same. The base time period was less 

in the case of predicted runoff. 



11000 13000 15000 17000 

Time (at 6min interval) 

Obs Pred Rain 


Figure 5.6; Best observed match in time series plot for 96 set 2 

The variation of RMSE with the increase in Kf values are presented in Figure 5.7. The 
RMSE is an absolute measure of the residual variance from the developed model 


Rain (mm) 
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Kf (mm/hr) 

Figure 5.7: Variation of RMSE with increase in Kf values 


'Fo decrease the peak values and to get rid of nunaber of peaks, the maximum infiltration 
rale w'as increased. Thirty trials were carried out by increasing the values from 10% to 
4000% with respect to infiltration values. For each case, RMSE and NSE were calculated. 
The variations of RMSE and NSE with the increase in maximum infiltration values are 
shown in Figure 5.8. It followed the expected behaviour to reduce the peak values by 
increasing maximum infiltration rates. After increment of 800%, peaks started 
disappearing. MBE was also calculated. The variation of MBE with the increase in 
Maximum infiltration values is shown in Figure 5.9. 




Mass Balance Error ^ ss RMSE 
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0 1000 2000 3000 4000 

Increament in Maximum Infiltration 


e 5.8: Variation in RMSE and NSC with increase in maximum infiltration 

values 



0 1000 2000 3000 4000 

Increament in Maximum Infiltration 


Figure 5.9: Variation in MBE with increase in Maximum infiltration values 


coei 
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5.4 Final 1996 trial 

The observed ground runoff were available for 10 months period (8.2.96 to 31.12.96). 
Using this data set, the runoff for 6 minute interval was generated. The runoff values were 
assumed to be uniform for the period in which data were missing. On the basis of previous 
discussion, the range of each variable was decided as follows: 

20 < Maximum infiltration < 40mm/hr 

0.5<Kf<50mm/lir 

3<Kst< 10 

Using these ranges, now a decision space has been generated. Sixty different trials were 
performed by moving points on each face of the cube. The time period for this study was 
from 00:00 on 13/8/96 to 12:30 on 31/12/96. But results obtained firom this study were not 
encouraging. So, the physics involved in each layer were studied to control the flow. 

5.4.1 Physics involved in soil moisture modelling 

The processes considered while developing WBM model needs to be understood to carry 
out further calibration trials. The entire processes are represented in Figure 5.10. The 
diagram is self explanatory. Deep percolation contributes to ground water which appears 
in channel, but the fate of deep percolation is not considered for the present study. 
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The significance of Kf and maximum infiltration values can be understood from Figure 
5.11. By defining Kf and maximum infiltration values as an input to model, the infiltration 
curve is developed by using well defined equations within the model. For a given soil 
moisture, the infiltration and percolation rates are calculated from the developed graph. At 
this point, it is important to note that this graph is developed only for the 1®^ layer. In 
present study, the soil is of 3 different types. So, out of 6 parameters used earlier to 
calibrate the model, only 3 are used now. It depends on which layer comes at the top. The 
model is instructed to use the saturated hydraulic conductivity and maximum infiltration 
values corresponding to the soil type, which exists at the top layer only. It is possible to 
shift the curve upward or downward by defining a factor. The Kf and maximum 
infiltration values will be multiplied by that factor and rest of the curve will be developed 
accordingly. 


( 



74 


Precipitation 

/ \ 

Intercepted precipitation Through fall precipitation 



Dropdown from leaves or stems of trees 




Figure 5.10: Processes considered in WBM model 




75 


Infiltration Rate 



Kr: Saturated hydraulic conductivity, WP: Welting Point, FC: Field capacity, MPV: Maximum pore volume 

Figure 5.11: Infiltration and percolation curve 

After realizing these facts, the following were decided to be checked for next calibration 
trials. 

i. Observe soil moisture profile besides time series graphs 

ii. Check the sensitivity of factor to be used 

iii. Significance of interflow should also be observed. 

Fifty different trials were carried out by changing the factor from 1 to 50. The factor of 4.0 
showed better results. For these trials, the Kf and maximum infiltration values were kept 
as mentioned below: 



16 


Kf (mm/hr) 

Maximum Infiltration 


(mm/hr) 

1.05 

1.15 

1.05 

1.15 

0.53 

0.58 


So, by defining a factor of 4.0, all the values will be multiplied by 4.0. Hence, it again 
showed values close to the best combinations decided since beginning (shown below) 
except the infiltration values for clay. 


Kf(mm/hr) 

Maximum Infiltration 


(mm/hr) 

4.5 

4.6 

4.5 

4.6 

0.05 

0.06 


It was observed that by increasing this factor, runoff increases. It was against the 
expectation that the runoff will decrease by increasing the Kf and maximum infiltration 
rate. One of the possible reasons could be back water effect. But since the ground water 
flow has not been considered for this study, so this reason does not seem logical. To 
include the effect of interflow, a channel of 5 m long was introduced artificially in the 
catchment at the last grid cell near to the outlet. This helps interflow to be added to the 


surface runoff. 
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5.4.2 Discussion on surface runoff and soil moisture profile 

It was decided to consider small duration for model run and investigate the output results 
in detail. The duration considered was from 12:00 on 26/8/96 to 23:54 on 30/9/96. The 
time series graph obtained is shown in Figure 5.12. The conclusions out of discussion in 
section 5.4.1 were fully incorporated during these model runs. It is clear from time series 
graph that except the presence of one high pealc, other peaks are almost matching for 
observed and predicted data. The problem of presence of such high peaks from predicted 
runoff values from model exists since beginning of calibration trials. So, at this juncture, it 
was felt that there could be error in observed ground data. It was urgent need to contact 
the concerned person who was directly involved during data collection. Dr. Andrew 
Western (the data were downloaded from his homepage) referred some publications which 
discuss such difficulty in developing model using available data for Tarrawarra catchment. 



Surface runoff (nim/h) 


78 



57239 58239 59239 60239 61239 62239 


Hme (at 6 min. interval) 

Obs Pred Rain 


Figure 5.12: Time series plot obtained from final trial 



1/4/96 29/4/96 27/5/96 24/6/96 22/7/96 19/8/96 16/9/96 14/10/96 


Figure 5.13: Time series plot obtained from literature 


Western [2000] mentions that the annual surface runoff during 1996 is overestimated by 
3 1 % from the model developed by his research group. The time series hydrograph showing 


Lam 
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these overestimations were collected from literature and is shown in Figure 5.13. By 
comparing the characteristics of hydrograph for the period 12:00 on 26/8/96 to 23:54 on 
30/9/96, it also shows presence of many peaks which are higher than observed ground 
data. 

i 

The next target was to analyse soil moisture pattern. Soil moisture pattern measured in 
ground on different dates were collected from literature and is presented in Figure 5.14. 
Western [2000] mentions that there are some small scale features that are critical to the 
runoff response of the catchment. These are the narrow bands of high soil moisture along 
the drainage lines and the model should be able to capture these. It also mentions that 
during winter, the catchment soil moisture is very high and extensive areas are saturated. 
The high moisture values in drainage lines on September 2"“^ and 20* are likely to be 
slightly too high due to surface water ponding introducing some measurement errors. The 
soil moisture profile obtained from WBM model at 12:00 on 2/9/96 is shown in Figure 
5.15. It is clear that the model is able to capture narrow bands of high soil moisture along 
the drainage lines. The remaining area with high value of soil moisture is due to 
approximations assumed during soil classifications in the catchment. 
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Figure 5.15: Soil moisture profile obtained from model 
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5.5 Calibration and validation using newly generated data set 

Based on the experience obtained at TU Darmstadt, another trial was made to achieve 
better calibration by trial and error in a more systematic way. The data generation and 
selection of calibration and validation periods are explained in section 4.8. The procedure 
to calibrate the model was kept same as discussed stepwise in section 5.1. The model 
calibration was tried for Period 4L, 4M and 4H. By using the best combination of 
parameters obtained during calibration, the model was validated for period 2L, 2M and 2H 
(described in section 4.8). The set of parameters obtained and corresponding time series 
graph are presented below from Figure 5.16a to Figure 5.16f. 



Figure 5.16a: Time series graph for calibration period of Low peak event i.e. 

period 4L 




Runoff (Cumecs) 
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Combination of parameters for period 4L 


Kf (mm/hr) 

Maximum Infiltration 

(mm/hr) 

Ks, 

6.5 

7.1 


6.5 

7.1 

0.3 

3.2 

3.5 




Time (at 6 min. interval) 


Figure 5.16b: Time series graph for validation period of Low peak event i.e. 

period 2L 
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0 
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Figure 5.16c: Time scries graph for calibration period of m edium peak event i.e. 

period 4M 


Combination of parameters for period 4M 

Maximum Infiltration Kjt 

(mm/hr) 

4.6 

4.6 3.0 

2.3 


Kf(mm/lir) 

4.5 

4.5 

2,25 
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Figure 5.16d: Time series graph for validation period of Medium peak event i.e. 

period 2M 
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Figure 5.16e: Time series graph for calibration period of High peak event i.e. 

period 4H 
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Combination of parameters for period 4H 

Kf(mm/hr) Maximum Infiltration Kst 

(mm/hr) 

2.5 2.8 

2.5 2.8 5.0 

1.3 1.4 



Figure 5.16f: Time series graph for validation period of High peak event i.e. 

period 2H 

5.6 Results and discussion 

The error statistics for calibration trials of different periods of year 1995 are presented in 
Table 5.2a. The date and time of these periods can be recalled from sections with same 
headings in this chapter. By looking at the errors, no definite pattern in any of the error 
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statistics are observed. Considering only R values suggests that 1995 trial - 1 with K^t of 
3.0 provides the best combination of parameters out of all trials. 

The error statistics obtained during calibration for data of different periods of year 1996 are 
presented in Table 5.2b. By looking at R values, it clearly indicates that the observed and 
predicted runoff are not matching. The AARE vaues are very high and other errors are 
also not considerable. However, only by visual interpretation of time series graph 
suggested best match in case of Final trial. 

The error statistics for calibration (Period - 4) and validation (period - 2) periods for low, 
medium and high peak flow events are presented in Table 5.2c. It is clear that the RMSE 
is decreasing from calibration to validation except in case of high peak flow event. 
Similarly R, AARE and SSE also showing improvement from calibration to validation for 
low and medium peak flow events. 

On the basis of time series graphs and error statistics obtained out of different calibration 
trials, reasonable results were obtained in case of 1995 trial - 1 with Kst of 3.0 and period 
4M, The combination of parameters in both the cases were averaged and the resulting 
parameter values are presented below. As discussed earlier in section 5.3.2, literature 
suggests that the Kf and Maximum infiltration values for clay should be half of the same 
for silt clay. 



Silt clay 
Silt 
Clay 
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Ki{mm/hr) 

Max. Inf. (mm/hr) 

K,t 

4.5 

4.6 


4.5 

4.6 

3.0 

1.175 

1.18 



As discussed earlier in section 5.3.2, literature suggests that the Kf and Maximum 

I 

infiltration values for clay should be half of the same for silt clay. By keeping this logic in 
mind the final best combination of parameters were declared as presented below. 


Kf(ram/hr) 

Silt clay 4.5 

Silt 4.5 

Clay 2.25 


Maximum Infiltration 
(mm/hr) 

4.6 

4.6 

2.3 


Kst 


3.0 


Table 5.2a: Error statistics of different calibration trials for 1995 


Period 

mSE 

NRM 

SE 

NSE 

El 

R 


SSE 

1995 

Trial -1, 
kst of 1.0 

0.000375 

2.973 


0.204 


87.13 

0.000065 

1995 

Trial- 1, 
kst of 3.0 




0.453 

0.945 

91.40 

0.000006 

1995 

Trial -1, 
kst of 5.0 


7.406 

-5.35 

-0.555 

0.869 

103.26 

0.000406 

1995 

Trial - 1, 
kst of 9.0 

0.000137 

10.896 

-12.75 

-1.374 

0.72 

123.58 

0.000878 

1995 

Trial-2 





0.0001 

135.8 

0.003169 































Tabic 5.2b: Error statistics of different calibration trials for 1996 


96 set 1 

0.002437 

Z706 

-0.045 

0.308 

0.247 

98.36 

0.000743 

96 set 2 

0.021042 

4.065 

-5.397 

-0.473 

0.122 

335.48 

0.12930 

Final trial 

0.014730 

3.710 

-TOlO 

-0.323 

0.437 

650.39 

0.08806 


Table 5.2c: Error statistics of different calibration and validation for 1996 


4L _ _ __ 

WEMSM 

0.794 


-0.099 

0.627 


0.003586 


0.000995 

1.3715 

0.396 

0.307 

0.686 

68.67 

0.00075 





0.1302 

0.719 

125.04 

0.01 


0.00375 

0.77 

0.6198 

0.301 

0.805 

83.175 

0.003502 


0.00786 

1.243 

0.595 

0.3099 

0.742 

215.54. 

0.01857 





-0.842 

0.0001 

496.32 i 

0.1745 


In the next chapter, the selected combination of parameter has been used to check the 
sensitivity of dynamic land use parameters. 






























Chapter 6 


Sensitivity analysis of dynamic land use parameters 

6.1 General 

In this chapter, first the different dynamic land use parameters have been briefly 
introduced. The seasonal variation of root depth, vegetation cover (BG), and leaf area 
index (LAI) were available for different type of forests. By applying the seasonal 
variations, the sensitivity of the parameters in terms of water balance have been observed. 

6.2 Dynamic land use 

Land use is based upon the function of land, i.e. the purpose for which the land is being 
used. A given land use may take place on one, or more than one pieces of land and several 
land uses may occur on the same piece of land. The word dynamic is added to stress upon 
the variations as a function of time e.g. seasonal cycle of vegetation and changes by 
introducing unproductive lands and vegetation of different types over the year. 

The catchment considered in the present study is covered of grass only and is used as 
pasteur land. The abovementioned land uses have been applied artificially. In the WBM, 
it is possible to alter the following parameters: 
i. Root depth 
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ii. Percentage of Vegetation cover (BG) and its seasonal cycle 

iii. Leaf area index (LAI) and its seasonal cycle 

iv. Unproductive land 
V. Roughness. 

The parameters are briefly described here. These informations are collected from Gurtz 
[2004]. 

The root depth 

The root depth is defined as the vertical distance form soil surface up to which 95% of all 
plant roots are found. Root growth is very sensitive to climate change and soil properties 
(especially moisture), so it is very important to properly quantity interaction between the 
climate, soil and plants. The soil moisture content available for evapotranspiration is 
controlled by the actual saturation of the soil and by the structure of the soil’s root zone. 

Leaf area index and vegetation density 

The seasonal cycle of vegetated surfaces can be well described through two seasonally 
varying parameters: The leaf area index (LAI) and the density of vegetation. In this way, 
it is possible to capture the different periodical behavior of the different classes of 
vegetation. LAI is the total area of leaves in relationship to the ground below them. 
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LAI describes a fundamental property of the plant canopy and its interaction with the 
atmosphere. Various soil — vegetation atmosphere models and hydrological models use 
LAI. The main hydrological processes related to seasonal development of these vegetation 
specific parameters are interception and evapotranspiration. Higher LAI means higher 
water storage capacity on leaves, which determine a higher potential for interception 
evaporation in the seasons with higher LAI. 

The surface roughness 

The average height of the different types of vegetation is assumed as index for the surface 
roughness. Roughness is used to calculate wind effects on evapotranspiration, which 
affects the rate of moisture transfer from the vegetation to the atmosphere. Hence, 
atmospheric and hydrological models require accurate estimates of surface roughness. By 
changing the values of these parameters, one can simulate different type of vegetation, 
their characteristics throughout the year and unproductive land. These parameters are 
linked to several processes of the hydrological cycle. These land uses alter the land in 
different ways. In terms of the movement of water in the basin, land use often changes the 
permeability of the land surface, thus affecting surface drainage systems and natural 
hydrology. Land use is important to the development of long term forecasts to guide 
present development. Systematic decision making regarding appropriate land uses for 
different parts of a watershed system has a central role in projecting or impacting valued 
landscapes. 



92 


6.3 Sensitivity of parameters for grass 

The effect of changes in root depth, BG, LAI and percentage coverage of unproductive 
land on water balance was observed for the period 00:00 on 30.03.96 to 23:54 on 30.12.96. 

6.3.1 Changes in root depth 

The root depth was varied from 0.45 mm to 0.8 mm. The graph showing the effect is 
presented in Figure 6.1. The precipitation is divided into through rainfall (TR), 
interception evaporation (IE), soil evaporation (SE), transpiration (T) and 
evapotranspiration (EPT). From Figure 6.1, it is clear that there is no effect in through 
rainfall, interception, and soil evaporation with the change in root depth. Transpiration and 
evapotranspiration values change by increasing root depth from 0.45 m to 0.7 m. In this 
case, BG was assumed to be 85 with LAI of 1.7 and Kst value was fixed as 3.0. 
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Root depth (m) 


Figure 6.1: Effect of changes in root depth 

The transpiration increases with an increase in the root depth up to 0.75 m. The maximum 
depth of soil considered in the present study is only 0.75 m. So, this might be the reason 
of constant transpiration value by increasing the root depth beyond maximum possible soil 
thickness. 


6.3.2 Changes in vegetation cover 


The change in water balance was observed for BG varying from 20 to 100. In this case, 
LAI was kept constant at 1,7. The effects are presented in Figure 6.2.. By increasing 
vegetation cover over the catchment, TR was decreasing because now a part of rain will be 
intercepted by leaves as expected. Hence, IE has increased. Since rain reaching ground 
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surface decreases with an increase in the BG; so SE is decreasing. The amount of T is 
increasing because of presence of large number of trees in the catchment. 



Figure 6.2: Effect of changes in BG 


6.3.3 Changes in leaf area index 

LAI were varied from 1.0 to 10.0. The effect on water balance is shown in Figure 6.3. 
TR, SE and EPT values are decreasing while IE and T values are increasing with the 
increase in LAI values. This looks very realistic that increment in LAI resulted into large 
obstruction in amount of rain reaching to ground and hence high IE. T is increasing due to 
presence of healthy or vegetated trees fully covered with leaves. 
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Leaf area Index 

Figure 6.3: Effect of changes in LAI 


6.3,4 Clianges in unproductive land 

For this study, the root depth, BG and LAI values were 0.45 mm, 85 and 1.7, respectively. 
The percentage of unproductive land in earlier cases was kept at 0. At this point, it is 
necessary to clearly define the meaning of percentage of unproductive land along with 
percentage of vegetation cover. In the WBM model, it is necessary to define ID 
corresponding to a part of the catchment. One can divide the catchment into two or more 
parts or it can be of single type only. By defining type, it means that the characteristics in 
terms of root depth, vegetation cover and LAI will be same for that portion of the 
catchment. These parameters are assigned to a single pixel and in case of single type of 
catchment, same characteristics is assigned to all the pixels of the catchment. 
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The partition within a single pixel can be well described using Figure 6.4. Initially the 
pixel is of unilorm characteristics. By defining a% of vegetation and b% of unproductive 
land, movement of water in a% of pixel will be slow due to presence of vegetation. The 
movement of water in b% of the rest (100-a) will be fast because in this portion, water is 
allowed to move horizontally not vertically. In the rest part of the pixel, the water 
movement will be normal. 



Figure 6.4: Partition of a pixel according to land use 

For the present study, the percentage unproductive land was varying firom 20% to 100% by 
keeping vegetation cover at 85% and single type to the whole catchment. The changes in 
water balance have been presented in Figure 6.5. The amount of IE, T and EPT has 
decreased by increasing percentage of unproductive land. Through rainfall is increasing, 
because obstruction by leaves is decreasing. By increasing percentage of improductive 
land within the catchment, automatically vegetated covered area is decreasing and hence in 
case of 1 00% unproductive land, there is no IE, SE, T and EPT. 
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0 20 40 60 80 100 120 

Unproductive land (%) 


Figure 6.5: Effect of introduction of unproductive land 

6,4 Different Land uses 

To observe the effect of different dynamic land use changes, 6 different types of forests 
were artificially generated in the catchment. Those are Young coniferous forest (NWl), 
Old coniferous forest (NW2), Young deciduous forest (LWl), Old deciduous forest 
(LW2), Young mixed forest (MWl), Old mixed forest (MW2). 

For each forest, two subcategories were generated to run the model. When BG is kept 
constant throughout the year and only seasonal variations of LAI is considered (denoted by 
‘a’) and when seasonal variations of both BG as well as LAI is considered (denoted by 
‘b’). 
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variations but there was no variation in BG as observed from Table 6.2. This resulted into 
small changes in transpiration values only. 


Table 6.2: Seasonal variations for coniferous forest 



NWla NWlb NW2a NW2b 

Types of coniferous tree 

Figure 6.6: Cumulative values of different components of rain for coniferous forest 

By comparing values from the graphs for young and old coniferous trees, it is cleat that 
TR, SE, T and EPT has decreased while IE has increased. This might have resulted due to 
high value of average LAI and BG in case of old coniferous tress. Time series graphs were 
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diawn for all the four cases (NWla, NWlb, NW2a, and NW2b). No definite conclusions 
could be diawn using these graphs. By manual observations, it was found that the 
variations were very small in all the components for different periods. 

6.4.1.2 Deciduous forest ^ 

The LAI and BG for deciduous forests have much variation with season. The values are 
presented in Table 6.3. In winter, it has LAI of 0 and in the summer seasons, it goes up to 
2.5. The BG for deciduous forest in winter is 10 for young trees and 15 for old trees, while 
in summer it reaches 60 for young trees and 70 for old trees. The cumulative values of 
different components of rain obtained from model are presented in Figure 6.7. 


Table 6.3: Seasonal variation for deciduous forest 
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LWla LWlb LW2a LW2b 


Types of deciduous tree 


Figure 6.7: Cumulative values of different components of rain for deciduous forest 


By comparing LWla and LWlb, it was observed that TR, SE and EPT has decreased while 
IE and T has increased. Average values of BG is kept constant as 32.5 for LWla. By 
introducing the seasonal variations of BG in case of LWlb, the decrement in through 
rainfall and EPT is obvious, because during the summer seasons, BG goes as high as 60 
and hence increment in IE. Similar pattern has been observed for LW2a and LW2b, which 
can be explained with similar logic. 


6,4. 1.3 Mixed forest 


The seasonal variations in LAI and BG for mixed forest are presented in Table 6.4 are 
nearly average of those for coniferous and deciduous forest. In winter season, the LAI for 
young and old trees is 1.4 and 2.0 respectively. In summer, for young trees it increases to 
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4.0 and for old trees it goes up to 5.54. The cumulative values of different components of 
rain for mixed forest are shown in Figure 6.10. 


Table 6.4; Seasonal variations for mixed forest 


Type 


Jan 

Feb 

Mar 

Apr 

May 

Jun 

Jul 

Aug 

Sept 

Oct 

Nov 

Dec 

Avg 

MWl 

BG 

42 

42 

42 

55 

65 

65 





42 

42 

52 







2.5 

4 

■ 

m 

DjU 

2.5 

1.7 

1.4 

2.5 

MW2 

BG 

46 

46 







65 


46 

46 

60 

MW2 

LAI 

2 




la 

5.54 

5.54 

5.54 

5.54 

4.2 


2.0 

3.6 



MWla MWlb MW2a MW2b 


Types of deciduous tree 

Figure 6,10: Cumulative values of different components of rain for mixed forest 

By comparing values for MWla and MWlb, it is clear that TR, SE, and EPT values are 
decreasing while IE and T are increasing. Similar is the case for MW2a and MW2b except 
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T has decreased in this case. From Table 6.4, it is clear that the BG is high in summer 
seasons. So, when variation in BG was considered, this resulted into increase in values of 
IE and T. 


/ 



Chapter 7 


Summary and conclusions 


7.1 Summary 

This thesis presents the findings of a study carried out with the purpose of using 
hydrologic models to check the sensitivity of dynamic land use parameters in a 
catchment. Many distributed hydrologic models were explored in an attempt to 
understand the use of RS and GIS in hydrological applications. 

The WBM model was used to carry out the distributed hydrologic modelling in this 
study. The climatic and hydrologic data of Tarrawarra catchment, Australia were used 
to apply WBM model. The available metereological data were processed to generate 
input files required to run the model. A method of trial and error was employed to 
estimate the seven parameters of the WBM model. The time series graphs between 
observed and predicted runoff and certain error statistics were used to check the 
accuracy of calibration. Based on various trials using data of different periods, a 
combination of parameters was selected for which there was good match in time series 
graph and errors were also encouraging. 

After fixing the parameters, the sensitivity of dynamic land use parameters were 
analysed. The seasonal change in root depth, percentage of vegetation cover in the 
catchment and leaf area index of the plant was considered as dynamic land use 



paiameters. Initially, the sensitivity of these parameters was checked for the catchment 
uniformly coveted with grass. In this case, the seasonal variations were not accounted. 
Latei, the change in water balance was observed by using seasonal variations of the 
land use parameters for coniferous, deciduous and mixed forest. 

7.2 Conclusions 

The results obtained during calibration of the WBM demonstrate that the WBM is able 
to capture the complex and non-linear relationship of rainfall runoff processes by 
taking into account the topography and spatial heterogeneity of soil type and land use. 
It is found that the model is sensitive to mainly seven parameters which change the 
model performance significantly. Moreover, the roughness is the most sensitive 
parameter. The non-uniform behaviour of infiltration rate with soil moisture played a 
major role in generating the response of the model. This created difficulty in 
estimating the values of parameters in the trial and error procedure. The time series 
graph and error statistics can only be able to demonstrate the direction to choose the 
values of parameters. It is unable to provide insights into the complex relationships. 
The findings of calibration presented in this report are only based on the evaluations of 
surface runoff. Due to discontinuity in available data, the soil moisture patterns could 
not be evaluated always except for a few occasions. 

The results of the sensitivity analysis of dynamic land use parameters show that the 
spatially distributed hydrologic model is a good tool to assess the hydrological change 
due to land use changes. GIS appears to be an efficient tool for presentation of 
seasonal variations of these parameters. Introduction of different types of vegetation in 



tlie catchment affects the water balance. This study clearly demonstrates that 
intei active integiation of spatial data and application of distributed hydrologic model in 
GIS environment provides a powerful tool for assessment of effects due to land use 
changes. Future changes in land use can also be incorporated in the niodel. Thus land 
use planning and management can be done efficiently. Based on the changes in 
different components of hydrologic cycle, it is concluded that the percentage change in 
water balance due to introduction of mixed forest lies almost at the middle of the 
changes observed due to coniferous forest and deciduous forest. 

7.2 Scope for future work 

No study is complete in itself, and there is always scope of its improvements. It is 
important to note that the primary analysis of the highly detailed data set is a critical 
stage in generating hydrologic models. In the present study, the periods used in the 
model were decided based on continuously available data set. It is required to fix the 
continuous 6 min. data set for at least one full year. This could be done by close 
observation of the pattern in the data set and accordingly fill the missing data. 

The evapotranspiration flux is one of the most significant fluxes in hydrologic 
modeling. In the present study, modified Penman-Monteith equation has been used to 
calculate evapotranspiration for 6 min interval. In most of the cases, the Penman- 
Monteith equation is used to find out daily evapotranspiration. So, the validity of the 
modified formula needs to be checked. 



At some occasion, the unavailable surface runoff values were generated by proper 
interpolation. For the surface runoff data measurement, water levels in the flume were 
measured every minute and were saved when the changes in water level exceed 1.5 
mm. This means that the missing data could be assumed as constant. Such methods 
may lead to errors in data that becomes limitation of the whole study. 

Use of optimisation algorithms to search the optimum values of saturated hydraulic 
conductivities can improve the model performance and needs to be investigated. In 
doing so, defining objective function is a major issue. By literature review and by 
gaining personal experience by assigning different error statistics as an objective 
function could lead to suitable objective function. Before defining constraints for 
optimisation problem, a detailed sensitivity analysis to find out possible optimal region 
is necessary. The studies presented in Duan (1992), Yapo (1998), Madson (2003) and 
Madson (2000) can provide necessary guideline in this direction. 

For the complete analysis of the effects of land use changes in the catchment, it is 
important to introduce urbanised area also to make it realistic urban development. This 
could be done by partitioning the catchment into different groups and then assigning 
different land uses for each group. This analysis will exactly simulate the 
establishment of buildings, trees and other land uses and its effect on future water 
availability inside and in surrounding area of the catchment. 

It would probably lead to improved calibration, if one were to calibrate the model 
against the soil moisture patterns and then use the model to predict runoff. Available 
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soil moisture patterns at key points in time could be used to estimate the field capacity 
and welting point. 
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Appendix - A 


Discussion on leaf area index 

LAI is a key parameter for global and regional models of biosphere/atmosphere exchange 
of carbon dioxide, water vapour and other materials. It also plays major role in determining 
the energy balance of the land surface. Methodologies for determining LAI have changed 
over the decades. It is useful to represent the inconsistencies that are actually found in the 
scientific literature. 

LAI is broadly defined as the amount of leaf area in a vegetation canopy per unit land area. 
LAI is a key structural characteristic of vegetation and land cover because of the role of 
green leaves in a wide range of biological and physical processes. LAI may be described 
most simply as 


LAI=^ 

G 

Where S is the functional (green) leaf area of the canopy standing on ground area G. 

Most commonly “S” is measured as the projected area. However, LAI ma y be more 
precisely defined in a number of different ways. For example, LAI may be measured as 
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the total surface area of leaves in a canopy. This will be equal to 2S for flat leaves and 
greater than 2S for needle shaped and succulent leaves and photosynthetic stems. 

According to Barclay (1998), there are at least four common measures of LAI. 

Def (1): Total LAI is based on the total outside area of the leaves, talcing leaf shape into 
account, per unit area of horizontal land below the canopy. 

Def (2): One-sided LAI is usually defined as half the total LAI, even if the two sides of 
the leaves are not symmetrical. 

Def (3): Horizontally projected LAI is the area of “shadow” that would be cast by each 
leaf in the canopy with a light source at infinite distance and perpendicular to it, summed 
up for all leaves in the canopy. 

Def (4): Inclined projected LAI represents the projected area of leaves taking into account 
individual leaf inclinations. 

Def (1) is rarely used. Def (2) suffers from the problem that the meaning of “one sided” is 
unclear for coniferous needles, highly clumped foliage, or rolled leaves.- Chen and Black 
(1992) suggest that the LAI of non-flat leaves should be defined as half the total 
intercepting area per unit ground area, and that def (3) should be abandoned. In the recent 
literature, def (4) is being used. 



Appendix B 


Soil horizons 


Most soils have a distinctive profile or sequence of horizontal layers. Generally, these 
horizons result from the soil processes of eluviation and organic activity. Five general 
layers are normally present in a typical soil: O, A, B, C, and R horizons (Figure Al). 


Eluviation 

Layer 

Illuviation 

Layer 



Figure Al: Typical layers found in a soil profile. 


The 0 horizon is the topmost layer of most soils. It is composed mainly of plant litter 
at various levels of decomposition and humus. Below it is the A-horizon. This layer is 
composed primarily of mineral particles, which has two characteristics: it is the layer in 
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which humus and other organic materials are mixed with mineral particles, and it is a 
zone of translocation from which eluviation has removed finer particles and soluble 
substances, both of which may be deposited at a lower layer. Thus the A horizon is 
dark in color and usually light in texture and porous. The A horizon is commonly 
differentiated into a darker upper horizon or organic accumulation, and a lower horizon 
showing loss of material by eluviation. 

The B horizon is a mineral soil layer which is dominated by illuviation. It receives 
material eluviated from the A horizon. This layer also has a higher bulk density than 
the A horizon due to its enrichment of clay particles. The B horizon may be colored by 
oxides of iron and aluminum or by calcium carbonate illuviated from die A horizon. 
These horizons can be divided into Al, A2 and A3 or Bl, B2 and B3 starting from top. 

The C horizon is composed of weathered parent material that has not been yet 
significantly affected by the pedogenic processes or translocation and organic 
modification. 

The R horizon consists of unweathered bedrock. 


Source: http://www.physicalgeography.net/fundamentals/10t.html 



